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The evolution equation for the flame surface 
density in turbulent premixed combustion 

By A. Trouve 1 AND T. Poinsot 2 


One central ingredient in flamelet models for turbulent premixed combustion is 
the flame surface density. This quantity conveys most of the effects of the turbulence 
on the rate of energy release and is obtained via a modeled transport equation, called 
the E-equation. Past theoretical work has produced a rigorous approach that leads 
to an exact, but unclosed, formulation for the turbulent E-equation (Section 1.2). 
In this exact E-equation, it appears that the dynamical properties of the flame 
surface density are determined by a single parameter, namely the turbulent flame 
stretch. Unfortunately, the flame surface density and the turbulent flame stretch 
are not available from experiments and, in the absence of experimental data, little 
is known on the validity of the closure assumptions used in current flamelet models. 
Direct Numerical Simulation (DNS) is the obvious, complementary approach to 
get basic information on these fundamental quantities. In the present work, three- 
dimensional DNS of premixed flames in isotropic turbulent flow is used to estimate 
the different terms appearing in the E-equation (Section 2.1). A new methodology is 
proposed to provide the source and sink terms for the flame surface density, resolved 
both temporally and spatially throughout the turbulent flame brush (Section 2.2). 
Using this methodology, the effects of the Lewis number on the rate of production 
of flame surface area are described in great detail and meaningful comparisons with 
flamelet models can be performed (Section 3). The analysis reveals in particular 
the tendency of the models to overpredict flame surface dissipation as well as their 
inability to reproduce variations due to thermo-diffusive phenomena. Thanks to the 
detailed information produced by a DNS-based analysis, this type of comparison 
not only underscores the shortcomings of current models but also suggests ways to 
improve them. 


1. Introduction 

1.1. The flamelet approach for turbulent premixed combustion 
Premixed turbulent combustion is the propagation of a chemical reaction zone 
through a turbulent, molecularly mixed region of fuel and oxidizer. The turbulent 
flame is characterized by the topology of the region in which reaction occurs: front, 
pockets, or large volumes. Depending on the relative values of various chemical 
and turbulence scales, dimensional analysis reveals a range of premixed combustion 
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inodes progressing from flamelets to distributed reaction zones to well-stirred reac- 
tors (Barrere 1974, Bray 1980, Borghi 1985, Peters 1986, Williams 1985, Poinsot et 
al. 1990). These modes correspond to different topologies of the reaction zone and 
require different approaches for both understanding and modeling. 

Experimental as well as theoretical evidence suggests that many technologically 
important flows occur in the flamelet burning mode. Flamelet combustion corre- 
sponds to chemical reaction occurring at fast time scales and short length scales 
relative to the turbulence. In this situation, the flame is confined to relatively thin 
layers within the turbulent flow field. 

In the flamelet regime, it is convenient to describe the flame-flow interactions in 
terms of two basic ingredients: a flame speed that characterizes the flame structure 
and the flame front surface area. For instance, the mean reaction rate may be 
written as the product of the mean fuel consumption rate per unit flame surface 
area times the flame surface density: 

{vr) = ( PuYr,u(Sc)s ) ( s '), ( 1 ) 

where u>r is the mass of fuel consumed per unit time and per unit volume; p u 
and Yr iU are respectively the density and the fuel mass fraction in the unburnt 
gas, Sc is the local integral of the reaction rate along the flame normal direction, 
Sc — f URdn, and characterizes the local combustion intensity; and E' is the flame 
surface area per unit volume. The flame surface density is defined as the expected 
value for E': E = (E'). 

In Eq.(l), the flamelet speed, { Sc)s , accounts for local variations of the reaction 
rate along the flame surface. Laminar flame theory indicates that the local flame 
structure is modified by flow divergence, usually characterized by the hydrodynamic 
strain rate acting in the flame tangent plane as well as by flame front curvature. 
Under certain conditions, these variations can become critical and lead to partial 
or total quenching of the flame. Recent studies, however, using Direct Numerical 
Simulations (DNS) suggest that quenching is a rather unlikely event for turbulent 
premixed flames (Poinsot et at. 1990). In addition, although the local combustion 
intensity may exhibit large variations along the turbulent flame front, particularly 
for non-unity Lewis number flames, DNS suggest that these variations always tend 
to cancel in the mean (Haworth & Poinsot 1992 1 _Rutland & Trouve 1991). In the 
simulations, the mean fuel consumption speed, Sc, defined as the area- weighted, 
space- averaged value of Sc integrated along the turbulent flame surface, remain 
within 10% to 30% from the one- dimensional, laminar flame speed value, si.\ 

Thus, it appears that in the absence of quenching, the mean fuel consumption 
speed, ~S C , is only weakly sensitive to the flow field and the principle effect of 
turbulence is for the fluctuating velocity field to wrinkle the flame and greatly 
increase its surface area. This phenomenon accounts for most of the increase in the 

f Note that Sc is a space- averaged quantity and should not be confused with ( Sc)s which is 
an area- weighted ensemble-average as defined in section 1.2 and, therefore, depends on location 
within the turbulent flame brush (see Figure 7) 
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overall burning rate due to the turbulence and is represented in Eq.(l) through the 
flame surface density, E. 

1.2. The evolution equation for the flame surface density, E 
In current flamelet models, the flame surface density, E, is obtained via a modeled 
transport equation. This equation was first postulated by Marble & Broad well 
(1977) based on phenomenological grounds. A more rigorous approach was later 
proposed by Pope (1988) and Candel k Poinsot (1990) who derive an exact balance 
equation for the flame surface-to-volume ratio, E': 


as; 

dt 


+ V.XE' = (V.X - nn : VX) E\ 


( 2 ) 


where X is the displacement speed of the flame surface, given by the sum of the fluid 
velocity and the flame propagation speed in the normal direction: X = u + urn; n is 
the unit vector normal to the flame surface; and where we use tensonal notations: 

(nn : VX) = 


riiTi 


dXj 
j dXj 


111 . V OX' , 

The right-hand side of Eq.(2) can also be expressed in terms of flame stretch. 
The flame stretch, fc, is defined as the rate of change of a Lagrangian flame surface 
element, 6 A: 

k - = QMl + X.V(6A) (3) 

dt dt 

A more useful expression for k is in terms of strain rate, flame curvature, and flame 
propagation speed (see for example Candel k Poinsot 1990): 


k = ar + 2 wk r 


(4) 


where a T is the rate of strain acting in the flame tangent plane: a T = V.u-nn : Vu; 
and k m is the flame surface curvature, as given by the divergence of the flame normal 
direction: 2k m = V.n. In Eq.(4), positive curvature is chosen convex towards the 

reset ant s 

Using Eq.(4), the balance equation for the flame surface-to-volume ratio can be 
re-written as: a . 

— + V.XE' = Jfc E' ( 5 ) 

dt 

When ensemble- averaged, this equation yields an exact balance equation for the 
flame surface density (Pope 1988, Cant et al. 1990): 


-jr- + V.(X) S E = {k)s E, (6) 

at 

where the flame surface mean of any quantity Q is given by: (Q)s = (QS')/( S') = 
(QE')/E. Note that surface means are different from standard means; in partic- 
ular, the surface mean of a quantity Q is different from the ensemble mean of Q 
conditioned on being at the flame location (see section 2.2). 
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Eq.(6) can be cast in various forms. For modeling purposes, it is useful to split the 
velocity vector into a mean component and a turbulent fluctuation: u = (U) + u'. 
We can then re-write Eq.(6) as follows: 

dE 

-I- V.(u ) S E + V.(wn) s E = { ar)s E + {At)s E + 2{wk m )s E, (7) 
where we use the following notations: 

( <*t)s = (V.u' - nn : Vu')s, 

( At)s = V(U) - (nn ) 5 : V(U) 

The three convective terms on the left-hand side of Eq.(7) are transport terms 
that correspond respectively to convection by the mean flow, turbulent diffusion, 
and flame propagation. The terms on the right-hand side of the equation are the 
source and sink terms for the flame surface density: {ot)s Is the turbulent strain 
rate acting in the flame tangent plane, (Ax)s is the strain rate due to the mean 
flow field, and 2 (wk m ) s is a term that accounts for the combined effects of flame 
curvature and flame propagation. 

The principle effect of turbulence is to increase the flame surface area, and (a T )s 
is without ambiguity a source term in the equation for E. The effect of the mean 
flow field as measured by (At)s is problem dependent; depending on the flow 
configuration, its sign can be positive or negative. We now focus attention on the 
last term in Eq.(7), referred to as the propagation term. 

In many situations, flame propagation effects merely counteract the wrinkling 
due to the turbulence and the propagation term, 2{wk m ) s , is, therefore, expected 
to be negative. Consequently, this term is usually described as a sink term in 
flamelet models. There are some situations, however, where this description is 
clearly incorrect. Since the propagation term includes some of the effects associated 
with intrinsic flame instabilities), this term must depend on the flame properties, 
thereby allowing for situations where its sign is positive and where the net effect 
corresponds to a production of flame surface. 

The exact importance of laminar flame instabilities for turbulent combustion is 
an open subject. Recent evidence, however, both experimental (Abdel-Gayed ct al. 
1984, Wu et al. 1990, Goix & Sheperd 1992) and numerical (Ashurst et al. 1987, 
Haworth &: Poinsot 1992, Rutland & Trouve 1991), suggests that the role of the 
Lewis number has been underestimated in the past. For instance, current flamelet 
models fail to account for the effects of the Lewis number on the rate of production 
of flame surface. The objective of the present study is to determine how present 
formulations might be improved to incorporate such effects. The approach is to 
analyze the source and sink terms in the equation for the flame surface density, 
with particular emphasis on how the Lewis number can affect their balance. This 
is accomplished using DNS, as described in the next section. 

f Using the terminology introduced to describe laminar flame instabilities, the propagation term 
represents the thermo-diffusive mechanism, while the strain term represents the hydrodynamic 
mechanism. These two instability mechanisms are coupled together and both account for Lewis 
number effects (Clavin 1985, Williams 1985) 
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2. Direct Numerical Simulation of turbulent premixed flames 
2.1 . Numerical method and configuration 
We use DNS to analyze the different terms appearing in the equation for the 
flame surface density. The simulations are performed using a three-dimensional, 
compressible Navier-Stokes solver that fully resolves the turbulent flow field. Spa- 
tial derivatives are computed with a modified Pade scheme that is sixth-order accu- 
rate (Lele 1990). Solutions are advanced in time using a third-order Runge-Kutta 
method (Wray 1990). Boundary conditions are specified using the NSCBC method 
(Poinsot & Lele 1992). Because of the otherwise prohibitive computational cost, 
simulations are limited to simple but finite-rate reaction schemes. In this work, the 
chemistry model is a single step, irreversible chemical reaction where the reaction 
rate depends exponentially on temperature (Arrhenius kinetics): 


wr = BpYr 



(S) 


where T a is the activation temperature and £ is a constant that depends on the 
flame speed. This formulation corresponds to a binary reaction in which one of the 
reactants, Yr, is strongly deficient as, for example, in fuel-lean combustion. Also, it 
is worth emphasizing that the simulations are not limited by the constant density 
assumption, and heat release effects are fully accounted for. 

Following Williams (1985), we re-write the reaction rate as: 

. . „ ( - 0(1 - 0 )\ ( 0 ) 

= \p\ R exp ^_ a( 1 _ 0 j ) ’ W 


where 0 is the reduced temperature, 0 = (T — T u )/{T b — T„); T u is the temperature 
of the fresh reactants; T b is the adiabatic flame temperature; and the coefficients 
A, a, and fl are, respectively, the reduced pre-exponential factor, the heat release 
factor, and the reduced activation energy: 

A = B exp(—/3/a), a = (T b — T u )/T b , and ft = aT a /T b (10) 


In the following, we use a = 0.75 and ft — 8. > 

Another important feature of the simulations is that transport coefficients are 
temperature dependent. These coefficients satisfy the following relations. 

p = p u (T/T u ) b , Le = X/pDc p = constant , Pr = pc p /X = constant , (11) 

where p, X, and D are the molecular diffusivities of, respectively, momentum, in- 
ternal energy, and species mass, b is a constant, and Le and Pr are respectively the 
Lewis number and the Prandtl number. We use b = 0.76, Pr = 0.75. Simulations 
have been performed for different Lewis numbers, Le = 0.8, 1.0, and 1.2. 

The selected computational configuration corresponds to a premixed flame em- 
bedded in three-dimensional, decaying, isotropic turbulent flow. The left- and right- 
hand sides of the computational domain are inflow and outflow boundaries while 
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periodic boundary conditions are applied at lateral walls. The calculations axe ini- 
tialized with fresh reactants on one side of the domain and burnt products on the 
other side; the two are separated by a plane laminar flame. Isotropic turbulence 
is initially located in the flow of fresh reactants, its velocity field being specified 
according to a model spectrum. The turbulence is characterized by a Kolmogorov 
length scale smaller than the thermal thickness of the laminar flame, ti^/St = 0.1 
where St = (T* - T u )/(dT/dx ) max ? and a turbulence intensity that is much higher 
than the laminar flame speed, u'/s L = 10. The initial turbulent Reynolds number, 
based on the Taylor microscale, is 50. The initial turbulent Reynolds number, based 
on the integral length scale, is 70. The grid resolution is 129 5 . 

The simulations describe the wrinkling of the flame zone due to turbulent motions 
as well as the combustion feedback due to dilatation and temperature- dependent 
transport properties. Note that the turbulence is decaying in time, and conditions 
are non-stationary. 

2.2 . Diagnostics 

All terms appearing in Eq.(7) may be obtained from the simulations. We now 
briefly describe how. The velocity vector and the velocity gradient tensor are readily 
obtained from the resolved flow field. To define flame-based quantities, we make use 
of concepts based on a thin flame picture. First, a progress variable, c, is introduced 
that is used to indicate location within the reaction zone, c = 1 — Yr, where Yr 
is the normalized fuel mass fraction. The progress variable varies monotonically 
through the flame from 0 in the reactants to 1 in the products. Constant progress 
variable surfaces may conveniently be used to define the flame front location: we 
use the surface c = Cf = 0.8. In addition, at any location on this surface, the local 
gradient of c defines the normal direction to the flame front: 


n = - 


Vc 
|Vc| ’ 


( 12 ) 


where n points into the fresh reactants. 

The propagation speed of the flame surface, w, is obtained from an expression 
analog to the well-known field equation (also called the G-ejjuation). Let us first 
consider a point on the flame surface, c — Cf. The velocity, X, at which this point 
must move to remain on the surface is given by: 

dc 

+ X.Vc = 0, (13) 

which, using Eq.(12), implies that: 


Y 1 

n “ 

and which yields the following expression for the flame propagation speed: 


w 


= X.n - u.n — 


1 Dc 
|Vc| DV 


(14) 


(15) 
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where quantities are estimated at the surface c = c f . Since in the simulations 
we solve for a conservation equation for the fuel mass fraction, Yr , , and since y 
definition c = 1 — Yr, the expression above can be readily obtained from the u 

The flame surface density, E, is a more subtle quantity. It includes both geomet- 
rical and statistical information. Following Pope (1990), the flame surface density 
is computed as the product of the expected value for the magnitude of the gradien 
of c, conditioned on being on the flame surface, times the probability of being on 

that surface: 


E' = |Vc| 6(c — Cf), and E = (S') = (|Vc| | c — c/) p(c/), (16) 


where p(c/) is the probability of c = c/. 

We now turn to the averaging problem. In the simulations, the flame brush 
propagates along the * direction, and the problem remains homogeneous m the 
y — z planes. Therefore, averaged quantities depend on x and time t only, and 
ensemble- averaging can be performed in the y — z planes. 

(Q){x,t) = yXr- f Q{x,y,z,t)dydz, (17) 

LyLt Z J 

where L y and L z are the y and 2 dimensions of the computation domain. The 
accuracy of this expression depends on the size of the computational domain with 
respect to the turbulent length scales. In the simulations, the integral length scale 
of the turbulent flow field grows as the kinetic energy decays; this growth, how- 
ever, is rather slow, and it was determined that the integral length scale remains 
at least 8 times smaller than L y and L z . Typically, in every y- z plane within the 
turbulent flame brush, the statistical sample consists of approximately 10 fully inde- 
pendent flame events, and, although we recognize that the statistics are somewhat 
undersampled, reasonable accuracy is expected when estimating the first moments. 

Conditional means are computed by integrating along the c = c/ contour: 

f Q dl 

(Q\c = c f )(x,t) = Jc ~ C/ d j~ ( 18 ) 

Jc—c/ 


Surface means are then obtained using the following relations. 


(Q)s(x,t) 


(QE') _ (QjVcj 1 c = cj) = L=c, 3l Vc l 
(E') (|Vc| \c = c f ) J c=cj |Vc| dl 


(19) 


Clearly, surface means differ from conditional means. 

The relations above provide a methodology to estimate the different means needed 
in our analysis. As shown in Eq.(16), the flame surface density, E, also requires 
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a) b ) 



Figure 1. A test case for flame stretch: the problem of a cylindrical, freely- 
expanding, premixed laminar flame: a) temperature contours; b) flame stretch: 
solid line based on Eq.(4); symbols based on Eq.(22). 


an estimate of the probability p(c/). Simple geometrical considerations lead to the 
following relation: 


p(cf) 


1 f dl 


( 20 ) 


and finally, Eqs.(16), (18) and (20) yield the following expression for E: 


; 1 [ £ |Vc| 

y(( If )> + (!;)> L,c, il 21 

Before applying these diagnostics to the turbulent flame simulations, we check 
the accuracy of our estimates for flame stretch and flame propagation speed using 
a model laminar flame problem as described in the next section. 


2.3. Validation of DNS-based estimates for flame stretch 

As seen in Eqs.(5) and (6), the flame stretch, A:, is the single relevant parameter 
that determines the growth rate of flame surface area. Since this growth rate is 
locally exponential, it is important to obtain accurate estimates for k . In our analy- 
sis, the flame stretch is obtained using Eq.(4). The overall accuracy of our analysis 
thus depends on our ability to predict correctly strain rate, flame curvature, and 
flame propagation speed. To check the accuracy of our estimates, we performed 
simulations of a cylindrical, premixed laminar flame expanding freely into an ini- 
tially quiescent medium (Figure la). In this model problem, the flame stretch can 
be directly measured from the growth of the flame radius: 


k 


1 dr / 

Tf dt 


( 22 ) 


where r f designates the radius of the flame contour, c = c/. Results based on this 
expression are compared with our DNS-based estimates in Figure lb. The very 
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FIGURE 2. DNS of turbulent premixed flames in isotropic flow. Snapshots of the 
flame surface after 3.5 turbulent eddy turn-over time: a) Le = 0.8; b) Le = 1.2. 
The flow is from top-left (reactants) to bottom-right (products). 


good agreement seen in Figure lb demonstrates that the flame stretch as well as 
the flame propagation speed can be accurately monitored with our diagnostics. 

3. Results and discussion 

As described in the previous section, the present study uses three-dimensional, 
direct numerical simulations of turbulent premixed flames in isotropic flow. Three 
different cases have been studied that correspond to turbulent flames characterized 
by the same laminar thermal thickness, 6t , the same laminar flame speed, si, 
embedded in the same initial turbulent flow field, but with different Lewis number, 
Le = 0.8, 1.0 and 1.2 (Figure 2). 

3.1. The overall effect of the Lewis number 
Figure 3 shows that the three cases exhibit large differences in the time history 
of the total reaction rate (space-averaged over the computational domain). After 4 
turbulent eddy turn-over time, t = 4r, the Le = 0.8 flame burns more than twice as 
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FIGURE 3. Lewis number effects on the overall combustion intensity. The total 
reaction rate is made non-dimensional by its initial value corresponding to a strain- 
free, plane laminar flame. Time is made non-dimensional by the initial, turbulent 
eddy turn-over time, t. 



Figure 4. Lewis number effects on the relative increase of total flame surface 
area. Time is made non-dimensional by the turbulent eddy turn-over time, r. 



Figure 5. Time evolution of the flame surface density, E, through the turbulent 
flame brush (reactants on the left; products on the right). Le = 0.8. E and x are 
made non-dimensional by the laminar thermal thickness, 6j \ Time is measured in 
units of the turbulent eddy turn-over time, t. 
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much as the Le = 1.2 flame. In agreement with the findings from previous studies 
(Ashurst et al. 1987, Haworth & Poinsot 1992, Rutland & Trouve 1991), the mean 
fuel consumption speed, Sc, irrespective of the value of the Lewis number, is found 
to be only a weak function of the turbulence: at all times, departures of Sc from 
the laminar value, s*,, remain within 20%. Thus, the principal effect of the Lewis 
number is to promote or inhibit the production of flame surface area. 

The drastic effect of the Lewis number on flame surface production is displayed 
in Figure 4 For Le = 1.0 and Le = 1.2, the flame surface area initially increases, 
reaches a maximum, and then decreases in time. The increase occurs as the turbu- 
lence wrinkles the initially flat flame surface. The flame then adapts to its turbulent 
environment, and, as the turbulence decays, the flame surface becomes smoother 
and relaxes to its initial state. The Le = 0.8 flame exhibits a strikingly different 
behavior: the flame surface area keeps increasing in time without saturation. Al- 
though saturation might be expected at later times, our simulations are limited by 
the size of the computational domain and this subsequent phase cannot be observed. 
In any case, the simulations indicate that saturation will not occur on a time scale 
characteristic of the turbulence, and, in that sense, the flame can be said to be 

unstable. „ . - A 

The differences between the Le = 0.8 and Le = 1.2 flames are in fact so pro- 
nounced that they can easily be observed by comparing instantaneous snapshots of 
the flame surface (Figure 2). For instance, for Le = 0.8, fingers of burnt products 
axe seen to propagate at a fast rate into the fresh reactants (Figure 2a). e e leve 
this ’’fingering” is an important ingredient of the flame instability process. Ihe 
’’fingering” is not observed in the Le — 1.0 or Le = 1.2 flames (Figure 2b). 

S.2. The source and sink terms tn the equation for E 
The effects of the Lewis number are now further studied by analyzing the structure 
of the terms appearing on the right-hand side of the equation for the flame surface 
density, E. As described in section 2.2, the analysis takes advantage of the fact that 
the problem is statistically one-dimensional and provides the source and sink terms 
for E as a function of time t and position x within the turbulent flame brush. 

Figure 5 compares several E-profiles through the turbulent flame brush taken at 
different instants in the simulations. The Lewis number is 0.8. At t = 0, E is a 
delta function located at x = 0. As time evolves, the turbulent flame brush gets 
thicker and propagates deeper into the reactants. Accordingly, the E-profile spreads 
out and shifts towards negative values of x. In the simulations, this shift is rather 
weak but can clearly be seen at the latest times ( t = 6r in Figure 5). Note that 
the integral of E through the flame brush gives the relative increase of total flame 

surface area: c / . . 


J E(x,f) dx = 


Sv(t) 


where Sv(t ) is the flame surface area within the computational domain of size V. 

The main advantage of the present analysis is to distinguish between the leading 
edge and the rear edge of the turbulent flame brush. The geometry as well as the 
dynamics of the flame differ quite significantly from one end of the reaction zone 
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Figure 6. Variations of the mean flame curvature, ( k m ) s , through the turbulent 
flame brush. Le = 0.8, t = 4 t. Flame curvature is made non-dimensional by the 
laminar thermal thickness, St . 



Figure 7 . Variations of the mean fuel consumption speed, (S c )s, through the 
turbulent flame brush. Le = 0.8, ( = 4r. S c is made non-dimensional by the 
laminar flame speed, si. 



FIGURE 8. Variations of the mean flame propagation speed, ( w)s , through the 
turbulent flame brush. Le = 0.8, t = 4 r. id is made non-dimensional by its value 
corresponding to a strain-free, plane laminar flame, w = sip(c = 0 )/p(c = c/). 
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to the other as shown in Figures 6 to 12. For convenience, in Figures 6 to 12, 
x-location within the turbulent flame brush is indicated by the local mean progress 

variable, (c). „ 

Figure 6 shows an instantaneous plot of the spatial variations of the mean flame 

curvature, { k m ) s , through the turbulent flame. The mean curvature goes from 
positive at the leading edge to negative at the rear edge of the turbulent flame. Since 
in non-unity Lewis number flames, the local combustion intensity, as measured by 
the reactant comsumption speed Sc, depends strongly on the flame curvature, and 
since the dependence is quasi-lineax (Haworth k Poinsot 1992, Rutland k Trouve 
1991), it might be inferred from Figure 6 that the statistical distribution of Sc is 
non-homogeneous as well. Figure 7 shows that this is indeed the case. For Le < 1, 

(. S c )s is a decreasing function of the mean progress variable, (c): the combustion 
intensity is higher at the leading edge than at the rear edge of the turbulent flame. 
For Le > 1, the trends are opposite: the combustion intensity is lower at the leading 
edge, close to (c) = 0, than at the rear edge, close to (c) = 1. For Le = 1, (S c )s 
remains approximately constant and equal to the laminar flame speed, si. 

In the flamelet regime, a flame element can be characterized by two speeds: Sc, 
which is a chemical rate, and w, which is a kinematic quantity and gives the velocity 
of the flame front with respect to the flow field (Eq.(15)). For a strain-free, plane 
laminar flame these two speeds are the same and equal tosf As pointed out by 
Poinsot et al. (1992), in the context of highly stretched flames, S c and w can be 
significantly different. Figure 8 shows the variations of the mean flame propagation 
speed, { w) S , through the turbulent flame. (w) s is an increasing function of the mean 
progress variable, (c). Comparison of Figures 7 and 8 indicate that, for Le = 0.8, 
the leading edge of the turbulent flame burns faster but propagates more slowly 
than the rear edge of the flame, which burns more slowly but propagates faster into 
the reactants 

It is worth emphasizing that the dynamical properties of the turbulent flame are 
not completely described by the knowledge of the distribution of the propagation 
speed, w, along the flame. To determine whether the flame surface will actually 
grow or contract, some information about the hydrodynamic flow field has to be 
included.! In other words, one needs to solve for the E-equation. 

We now turn to the terms appearing on the right-hand side of the equation 
for E. While the strain term, ( a T )s , remains approximately constant through 
the turbulent flame (Figure 9), the propagation term, 2(ti>*m)s, exhibits strong 
variations and decreases from positive values on the unburnt side, close to (c) — , 

to negative values on the burnt side, close to ( c ) = 1 (Figure 10). The net effect on 
the surface growth rate is given by the flame stretch, ( k)s — ( a T)s + 2 (w ’m)s- 
Figures 9 and 10 indicate that both contributions to stretch have the same order of 

magnitude. 


f This is best seen in the context or laminar flame instabilities, where the classical linear theory 
shows that the stability problem is not solved at the level of determining the Markstein length but 
also requires solving for a dispersion relation, which includes hydrodynamic effects (Clavin 1985, 
Williams 1985) 
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FIGURE 9. Variations of the mean strain rate, (clt)sj through the turbulent flame 
brush. Le = 0.8, t = 4 r. Strain rate is made non-dimensional by the characteristic 
flame time, 8t/sl- 



Figure 10. Variations of the mean propagation term, 2{wk m )s , through the 
turbulent flame brush. Le = 0.8, t = 4 r, wk m is made non-dimensional by the 
characteristic flame time, 6t/sl- 



FIGURE 11. Variations of the mean flame stretch, (k )s, through the turbulent 
flame brush. Le = 0.8, t = 4 r. A: is made non-dimensional by the characteristic 
flame time, 6t/$l ■ 
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Figure 12. Lewis number effects on the variations of the mean flame stretch, 
(1c) s , through the turbulent flame brush, and comparison with flamelet models. 
t = 4r. k is made non-dimensional by the laminar flame time, S T /s L . 


Figure 11 presents the mean flame stretch profile through the turbulent reaction 
zone. At the leading edge, strain rate and flame propagation effects are locally 
cumulative and the overall balance is strongly positive. The leading edge of the 
turbulent flame is a region of strong production of flame surface area. On the 
contrary, the propagation term takes large negative values on the burnt side. In 
that region, strain rate and flame propagation effects are locally opposite, and the 
overall balance is negative. The rear edge of the turbulent flame thus appears as a 
region where flame surface area gets strongly dissipated. ... 

Figure 11 spatially resolves the balance between production and dissipation of 
flame surface area. The net effect is given by defining a mean stretch, k, space- 
averaged throughout the flame brush: 


j _ df Svjt) , 

dt y L y L z ’ 



(24) 


k gives the instantaneous rate of change of the flame surface area in the compu- 
tational domain. If k is positive, the flame surface grows; if negative, the flame 
surface contracts. The next section further discusses the effect of the Lewis number 
on the spatially-resolved flame stretch profile, as well as the resulting impact on the 
net mean flame stretch, fc, and presents some comparison with flamelet models. 

3.S. Comparison of DNS results with flamelet models 

Figure 12 compares the mean flame stretch profiles, ( k)s , plotted for different 
Lewis numbers. In all cases, stretch takes large negative values on the burnt side, 
close to (c) = 1. The effect of the Lewis number is not visible in that region. On 
the contrary, at the leading edge of the flame, there are large differences between 



340 


A. Trouve & T. Poinsot 


the different Lewis number cases. This suggests that the turbulent flame is most 
sensitive to Lewis number effects on the unburnt side, whereas it remains unaffected 
on the burnt side. Using Eq. (24), we can quantify the overall differences previously 
observed in Figure 4; we find (in units of the laminar flame time): k = 2.6 for 
Le = 0.8, and the flame surface area is quickly growing; k ~ 0.0 for Le = 1.0 and 
Le = 1.2, and the flame surface area remains roughly constant. 

Also plotted in Figure 12 is a comparison with flamelet models. Many current 
flamelet models use a transport equation for the flame surface density. Different for- 
mulations of this equation have been proposed in the literature (Marble &c Broad well 
1977, Candel ct al. 1990, Cant et al. 1990, Borghi 1990). For the sake of compari- 
son, we use a closure assumption similar to the one proposed in the Coherent Flame 
Model by Marble & Broadwell (1977) and Candel et al. (1990). In this formulation, 
the turbulent flame stretch is written as: 

< 25 > 

where k t is the turbulent kinetic energy and e its dissipation; (V/*) is the ensemble- 
averaged fuel mass fraction. The first term on the right-hand side of Eq.(25) rep- 
resents straining due to the flow motions and is assumed to scale with the integral 
time scale of the turbulence; the second term is a disparition term associated with 
flame propagation and is assumed to scale with the laminar flame speed, s Ly and 
the flame surface density, E. 

Figure 12 shows that this model is indeed able to reproduce qualitatively the 
spatial structure of the balance between production and dissipation of E, going 
from production at the leading edge of the turbulent flame to dissipation at the 
rear edge. However, the 1 /(Yr) behavior of the disparition term leads to numerical 
difficulties on the burnt side of the flame. The model, therefore, overpredicts the 
dissipation of E, near (c) = 1, and gives a negative mean flame stretch, k = -3.6, 
in strong disagreement with the values reported above. In addition, the disparition 
term in Eq.(25) is always and everywhere negative and cannot account for the 
possible transition to unstable flame conditions as observed in the simulations. 


4. Conclusion 

Flamelet models constitute one of the most common approach for turbulent pre- 
mixed combustion. In these models, the flame surface density is a central ingredient 
that conveys most of the effects of the turbulence on the rate of energy release. The 
flame surface density is usually obtained via a modeled transport equation, called 
the E-equation, first postulated by Marble k Broadwell (1977). Recent theoretical 
work, based on conservation equations for surfaces and volumes in a turbulent flow 
field, has produced a more rigorous approach that leads to an exact, but unclosed, 
formulation for the turbulent E-equation (Pope 1988, Candel k Poinsot 1990). In 
this exact E-equation, it appears that the dynamical properties of the flame sur- 
face density are determined by a single parameter, namely the turbulent flame 
stretch. Unfortunately, the flame surface density and the turbulent flame stretch 
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are extremely difficult to measure and are simply not available from experiments. 
Therefore, little is known on the validity of the closure assumptions used in current 
formulations of flamelet models. 

Direct Numerical Simulation (DNS) is the obvious, complementary approach to 
get basic information on these fundamental quantities. In the present work, three- 
dimensional DNS of premixed flames propagating in isotropic turbulent flow have 
been used to estimate the different terms appearing in the E-equation. A new 
methodology has been proposed to provide the source and sink terms for the flame 
surface density, estimated as a function of time and position within the turbulent 
flame brush. Using this methodology, the effects of the Lewis number on the rate 
of production of flame surface area are described in great detail. Principal findings 
axe that - (1) the balance between production and dissipation of flame surface area 
is strongly non-homogeneous: the leading edge of the turbulent flame is a region 
of production of flame surface area, whereas the rear edge is a region where flame 
surface gets strongly dissipated; (2) the turbulent flame is most sensitive to Lewis 
number effects at the leading edge, whereas it remains unaffected on the burnt side. 
These results suggest that most of the important dynamical features of turbulent 
flames take place at the leading edge of the reaction zone. 

Detailed comparisons with flamelet models were also performed. The analysis re- 
veals the tendency of the models to overpredict flame surface dissipation as well as 
their inability to reproduce variations due to thermo-diffusive phenomena. Thanks 
to the detailed information produced by a DNS-based analysis, this type of compar- 
ison not only underscores the shortcomings of current models but also suggests ways 
to improve them. Future work will focus on the development of a new formulation 
of the E-equation that would incorporate thermo-diffusive mechanisms. 
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